#ifndef _BRESENHAM_H_
#define _BRESENHAM_H_

#include <cmath>


class Bresenham {

public:

  template <typename Func>
  static void line2D(int x0, int y0, int x1, int y1, Func &func) {
    /* Swap x and y if the absolute value of the slope
     * is greater than one. */
    if (std::abs(y1 - y0) <= std::abs(x1 - x0)) {
      line2D(x0, y0, x1, y1, x0, y0, func);
    } else {
      line2D(y0, x0, y1, x1, x0, y0, func);
    }
  }

  template <typename Func>
  static void line3D(int x0, int y0, int z0, int x1, int y1, int z1, Func &func) {
    int dx = std::abs(x1 - x0);
    int dy = std::abs(y1 - y0);
    int dz = std::abs(z1 - z0);

    /* Check for the direction with the biggest extension. */
    if ((dx >= dy) && (dx >= dz)) {
      line3D(x0, y0, z0, x1, y1, z1, x0, y0, z0, func);
    } else if ((dy >= dx) && (dy >= dz)) {
      line3D(y0, x0, z0, y1, x1, z1, x0, y0, z0, func);
    } else {
      line3D(z0, x0, y0, z1, x1, y1, x0, y0, z0, func);
    }
  }

private:

  template <typename Func>
  static void line2D(
      int &x0, int &y0,
      int x1, int y1,
      int &rx, int &ry,
      Func &func) {
    int xInc = (x1 > x0) ? +1 : -1;
    int yInc = (y1 > y0) ? +1 : -1;

    int dx = std::abs(x1 - x0);
    int dy = std::abs(y1 - y0);
    int d = 2 * dy - dx;
    int incE = 2 * dy;
    int incNE = 2 * (dy - dx);

    while (func(rx, ry) && (x0 != x1)) {
      if (d <= 0) {
        d += incE;
      } else {
        d += incNE;
        y0 += yInc;
      }
      x0 += xInc;
    }
  }

  template <typename Func>
  static void line3D(
      int &x0, int &y0, int &z0,
      int x1, int y1, int z1,
      int &rx, int &ry, int &rz,
      Func &func) {
    int xInc = (x1 > x0) ? +1 : -1;
    int yInc = (y1 > y0) ? +1 : -1;
    int zInc = (z1 > z0) ? +1 : -1;

    int dx = std::abs(x1 - x0);
    int dy = std::abs(y1 - y0);
    int dz = std::abs(z1 - z0);
    int d1 = 2 * dy - dx;
    int d2 = 2 * dz - dx;
    int incE1 = 2 * dy;
    int incNE1 = 2 * (dy - dx);
    int incE2 = 2 * dz;
    int incNE2 = 2 * (dz - dx);

    while (func(rx, ry, rz) && (x0 != x1)) {
      if (d1 <= 0) {
        d1 += incE1;
      } else {
        y0 += yInc;
        d1 += incNE1;
      }
      if (d2 <= 0) {
        d2 += incE2;
      } else {
        z0 += zInc;
        d2 += incNE2;
      }
      x0 += xInc;
    }
  }
};

#endif /* _BRESENHAM_H_ */

